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ABSTRACT 

We present the first radiation-hydrodynamic model of a protoplanetary disc irradi- 
ated with an X-EUV spectrum. In a model where the total ionizing luminosity is 
divided equally between X-ray and EUV luminosity, we find a photoevaporation rate 
of 1.4xlO~* M0yr^^, which is two orders of magnitude greater than the case of EUV 
photoevaporation alone. Thus it is clear that the X-rays are the dominant driving 
mechanism for photoevaporation. This can be understood inasmuch as X-rays are ca- 
pable of penetrating much larger columns lO^^cm"^) and can thus effect heating 
in denser regions and at larger radius than the EUV can. The radial extent of the 
launching region of the X-ray heated wind is 1-70AU compared with the pure EUV 
case where the launch region is concentrated around a few AU. When we couple our 
wind mass-loss rates with models for the disc's viscous evolution, we find that, as 
in the pure EUV case, there is a photoevaporative switch, such that an inner hole 
develops at ~ 1 AU at the point that the accretion rate in the disc drops below the 
wind mass loss rate. At this point, the remaining disc material is quickly removed in 
the final 15-20% of the disc's lifetime. This is consistent with the lO^yr transitional 
timescale estimated from observations of T-Tauri stars. We however note several key 
difference to previous EUV driven photoevaporation models. The two orders of mag- 
nitude higher photoevaporation rate is now consistent with the average accretion rate 
observed in young stars and will cut the disc off in its prime. Moreover, the extended 
mass-loss profile subjects the disc to a significant period (^20% of the disc's lifetime) 
of 'photoevaporation starved accretion'. We also caution that although our mass-loss 
rates are high compared to some accretion rates observed in young stars, our model 
has a rather large X-ray luminosity of 2 x lO'^^erg s~^; further modeling is required 
in order to investigate the evolutionary implications of the large observed spread of 
X-ray luminosities in T-Tauri stars. 

Key words: accretion, accretion discs - circumstellar matter - planetary systems: pro- 
toplanetary discs - stars: pre-main-sequence - X-rays: stars. 



1 INTRODUCTION 

The structure and evolution of protoplanetary discs is cur- 
rently an important and much debated topic in astrophysics. 
Detailed knowledge of disc evolution is needed to constrain 
the conditions under which planets form. In particular, an 
important issue in disc evolution is the mechanism by which 
the disc is eventually dispersed and the timescales involved, 
which consequently sets the time available for planet forma- 
tion to occur. 

For solar type stars it is now well established obser- 



vationally that at an age of a few Myr a large percent- 
age of stars are still surrounded by discs that are op tically 
thick at infrared wavelengths (e.g.^Ha isch et aLll200ll ). How- 
ever after ~ lO^yr most stars are observed to no longer be 
surrounded by discs. The relatively small number of discs 
that are observed to be in transition (between being sur- 
rounded by an optically thick disc and having no detectable 
disc) implies that the transition occurs on a much shorter 
~ lO ^yr timescale l|Kenvon fc Hartmann|[l995l : iDuvert et al.l 
l200d ). This 'two timescale' phenomenology provides a strong 
constraint on models of disc evolution and dispersal and is 
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inconsistent with standard viscous-accretion mode ls which 
predict a power-law decline (jHartmann et allllQQSl 'Fl. 

One of the most successful models in reproducing this 
'two-timescale' behavior is photoevaporation by radiation 
from the central star. Recent detection of [ Ne ll] emis- 
sion l ines (iPascucci 2OO7I: Herczeg et al.l 200?!: iNaiita et alj 



I2OO9I : iFlaccomio et al.ll2009l : iPascucci et ahlbOOgf T has been 
interpreted as evidence for irradiation of the disc's at- 
mosphere by energetic photons (extreme ultrav iolet or X- 
ray) and possibly associate d photoevaporation (jAlexandeil 
l2008bl : |P ascucci et al.ll2009l l. However there is still no con- 
sensus over whether photoevaporation is best explained in 
terms of far ultraviolet (FUV, hv < 13.6eV), extreme ul- 
traviolet (EUV, hu > 13.6eV) or X-ray irradiation (hi' > 
lOOeV), or some combination of the above. EUV photoevap- 
oratio n (of optically thick discs) is theoreticall y well devel- 
oped l|Hollenbach et al.lll994l : iFont et al.l[2004l ) with a pre- 
dicted mass-loss rate of ~lO~'^°M0yr~^. The mass-loss pro- 
file is dominated by a region within several AU of the char- 
acteristic (gravitational) radius rg = 8.9(M«/Mq)AU where 
the internal energy of the gas exceeds the escape potential. 
Co upling this mass-loss profile to viscous evolution mod- 
els ([Clarke et al...2001i ). this 'EUV-switch' model shows that 
once the accretion rate through the disc drops to the pho- 
toevaporative mass-loss rate, a hole develops and the inner 
disc drains onto the central star on the short viscous accre- 
tion timescale at a few AU. The inner region is now optically 
thin to EUV photons and a highe r 'direct- wind' mass-los s 
rate of ~lO"^M0yr"^ is obtained l|Alexander et al.ll2006al ). 
The disc then phot oevaporates from the i nside out on a 
timescale of ~10^yr (I Alexander et al]|2006lj ). The relevance 
of the EUV-switch model, however, depends, albeit weakly, 
on the EUV luminosity of the star; the canonical mass loss 
rates quoted above correspond to an EUV luminosity of 10''^ 
photons s~^. This cannot be verified by direct observation 
owing to the strong absorption of ionizing photons by in - 
tervening neutral hydrogen. While [Alexander et al.| (|2005| ) 
derived EUV fiuxes from young solar type stars which are 
consistent with those required by photoevaporation, it has 
been pointed out (e.g. Glas sgold et al.ll2004l : lErcolano et al.l 
[2009b; Glassgold ct al. 200^) that the EUV photons emitted 
by the star do not necessarily reach the disc, given that EUV 
photons are extremely easily absorbed, by small columns of 
neutral material surrounding the young star. X-ray obser- 
vations of young T-Tauri stars show that neutral hydrogen 
columns of up to lO^^cm"'^ are common around many ob- 
jects, probably due to accretion funnels, implying that EUV 
photoevaporation is unlikely to be effective until late times. 

FUV driven neutral fiows have also been considered as 
a possible mechanism for dis c dispersal. In rich cluster en- 
vironments, Ijohnstone et al. I ll99S') showed that FUV ra- 
diation from nearby stars can dominate the photoevapora- 
tion rate. This, however doesn't explain the rapid disper- 
sal of dis cs around young sta r s with no nearby source of 
FUV fiux. iGorti fc HoUenbachl (|2009^ have considered FUV 
driven fiows from the central star and argue that mass-loss 



It is a matter of current debate wheth er this two-timescale 
description also applies to low mass stars (see|SiciIia-Ae:uiIar et al] 



rates of order 10 M0yr can be obtained at large radius 
(> lOOAU). 

In contrast to EUV radiation, X-rays are able to 
penetrate larger columns of neutral material. Also X-rays 
are largely unaffected by accretion funnels (since they 
are emmitted in the magnetosphere, above the accretion 
funnels); furthermore the X-ray lu minosity function of 
young T-Tauri star is well k nown l|Preibisch et al.l l2005l : 
lAlbacete Colombo et al.l l2007'). For these reasons several au- 
thors have considered X-ray irradiati on as a possible drivini 
mechanism for photoevapo r ation (e.g. Alexander et al.l l200' 
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description also applies to low mass stars lsee|bi( 
I2OO8I: lErcolano et alll2009al : ICurrie et aUbOOgli ^ 



iGorti fc Ho Uenbach' '200g|; 

While lAlexandcr ct al.^ l|2004 ) and _ 

(|2009l ) quantitatively disagree on the mass-loss rates due 
to X-ray irradiation, they both conclude that X-ray heat- 
ing is conside r ably l e ss imp ortant than EUV heating. 
lErcolano et~al] (|2008bl . l2009bh obtained mass-loss rates in 
the range 10~* — lO~^''M0yr~^ and argue that the X-ray 
heating will dominate the mass-loss and hence the evolution 
of the disc. While the hydrodynamics of EUV photoevapo- 
ration have been solved in detail, there ar e no prior hydro- 
dynamical studios of X-ray irradiated gas. [Alexander et al.l 
(2004); Gorti & HoUonbach (2009,) both_useal+lD descrip- 
tion o f the radiative transfer, and [Ercolano et al.[ (|2008bl. 
[2009bl ) uses a 3D approach; however, in all previous studies 
the disc structures are assumed to be in hydrostatic equilib- 
rium and the mass loss-rates are obtained either by assuming 
a sonic fiow from the location where the local gas temper- 
ature equals the escape temperature (the pCs meth od, e.g. 
[Alexander et all [2004 [Ercolano et aLli2008bl . [2"009bD . or by 
coupling a spherical Parker wind solu tion to the hydrostatic 
structure (jCorti fc Holle nbach 2009). While these methods 
may (or may not) give a reasonable order of magnitude es- 
timate of the total mass-loss rate, they are unlikely to give 
a reliable radial profile of the mass-loss rate, as is needed 
to determine the evolution of viscously evolving photoevap- 
orating discs. 

In this paper we c ouple the radiative transfer models of 
[Ercolano et al.[ l|2009bl ) to two-dimensional hydrodynamical 
models to obtain a self consistent fiow solution for a primor- 
dial disc (i.e. an optically thick disc extending to the dust 
destruction radius) and for inner hole discs (i.e. gas discs 
with cleared inner holes). In section 2 we present our ap- 
proach to the radiation-hydrodynamics. In sections 3 and 4 
we apply these methods to a primordial disc and discs with 
inner holes and present the basic results. In section 5 we 
show the results of various tests we conducted. In section 
6 we present the results of coupling our wind profiles to a 
viscous evolution model. In section 7 we compare our work 
to previous results and discuss the model's limitations in 
section 8. In section 9 we discuss the observational conse- 
quences of our model and summarise our main findings in 
section 10. 



2 METHODS 

In what follows we undertake 2D hydrodynamic calculations 
in which - under the assumption that the gas is in radiative 
equilibrium - we parametrize the gas temperature as a func- 
tion of local variables. This parametrization is based on the 
results of 3D Monte Carlo radiative transfer calculations us- 
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ing MOCASSIN l|Ercolano et al.ll2003l , |2005| . l2008al ). We then 
verify a posteriori that the cooling timescale is less than the 
flow timescale and hence justify our assumption of radiative 
equilibrium. Below we set out the details of the radiative 
and hydrodynamic elements of the calculation. 

2.1 Radiative Transfer 

lErcolano et~al] l|2008bl . l2009lj ) employed the three- 
dimensional Monte Carlo photoionization and dust 
radia t ive tra nsfer code, MOCASSIN (Ercolan o et al.l |2003| . 
I2OO5I . l2008al ). to calculate the temperature and ioniza- 
tion structure of a typical T-Tauri disc irradiated by 
chromospheric X-ray and EUV flux. MOCASSIN uses a 
stochastic approach to the transfer of radiation and allows 
for arbitrary geometries and density distributions. The 
version of the c ode used for t hese models was modified by 
lErcolano et al. I (|2008bl , [2009bh and we refer to these articles 
for a detailed description of the co de and model se t up. We 
use the model labeled FS0H2Lxl by lErcolano et all l|2009bl ') 
and summarise here only the basic input parameters. 

2.1.1 Intial Disc Configuration 

The initial density structure of a protoplanetary disc sur- 
rounding a 0.7 Mq star (Teff = 4000 K, R, = 2.5 R©) was 
taken from the iD'Alessio et al.l (1200 ll ) set and was chosen 
as the model that best fit the median SED in Taurus. This 
model is in vertical hydrostatic equilibrium in the case that 
the dust temperatures are set by irradiation by the optical 
radiation from the central star and under the assumption of 
full thermal coupling between the gas and dust. The total 
mass in the disc is 0.026 M© with an outer radius of 500 AU. 

This disc structure, which provides our starting density 
distribution, employs a bimodal dust distribution, where at- 
mospheric dust follows the standard MRN model but with a 
very low dust to gas mass ratio (2.5xl0~^ for graphite and 
4.0x10"^ for astronomical silicates) and interior dust con- 
sists of larger grains with a size distribution still described 
by a power law of index -3.5, but with Omm = 0.005 /im and 
Omax = 1 mm and a dust-to-gas mass ratio of 3.1 x 10~^ 
for graphite and 5.0 x 10~^ for astronomical silicates. The 
transition between atmospheric and interior dust occurs at 
a height of 0.1 times the midplane gas scale height. 

High-energy dust absorption and scattering coefficients 
are calculated fr om the dielect ri c con stants for graphite 
and silicates of iLaor fc Draind l| 19931 ). which extend to 
the X-ray domain. Spherical grains are assumed and we 
use the standard Mie scattering series expansion for com- 
plex refractive indices, x\m\ < 1000, where x = 2a/ X 
is the scattering parameter (see iLaor Draind [l993l ). For 
x\m\ > 1000 and x\m — 1| < 0.0 01 we use Rayleigh-Gans 
theory (|Bohren fc Huffmanl \l983i ). and for x\m\ > 1000 
and x\m — 1| > 0.00 1 we use the treatment specified by 
iLaor fc Draind (|l993l ). which is based on geometric optics. 

2.1.2 Method 

This initial disc was then irradiated by a synthetic coronal 
spectrum of X-ray luminosity Lx = 2x lO^^erg s~^ (see Fig- 
ure [T] ) and total ionizing luminosity of Ltot = 4xl0^''erg 
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Fig ure 1. Input sp e ctrum used in MOCASSIN calculation. See text 
and lErcolano et al.l l|2009bl ) for details. 

s~^. This synthetic spectrum wa s a thermal spectr u m gen - 
erated by the plasma code of iKashvap fc Drakd (|2000h : 
the emission measure distributi on is based on that d erived 
for RS CVn type binaries by ISanz- Forcada et al] (l2002i ) 
(which peaks a t 10'*K) and fi ts to Chandra spectra of T 
Tauri stars by iMaggio et al.l (|2007h ( which latter peaks 
at around logT = 7.5K). Our use of the input spectrum 
shown in Figure [T] implies that we assume a zero column 
of neutral material close to the s ource. The experiments de- 
scribed in lErcolano et"al] l|2009bl ) however show that screen- 
ing columns of up to 10^" cm~'^, which would only absorb 
the EUV component of the spectrum, have little effect on the 
resulting photoevaporation estimates; we thus expect our re- 
sults to generalise to the case of moderate screening (< 10^" 
cm~'^) since this has little effect on the X-rays in the keV 
range whi ch dominate the disc h eating (see below). As de- 
scribed in lErcolano et af] (|2009bl ). the disc temperatures are 
updated in response to irradiation by the X-EUV spectrum 
and the density structure updated so as to attain a situation 
of hydrostatic equilibrium for this updated temperature pro- 
file. These temperature and density distributions are then 
updated iteratively until a situation of both thermal and 
hydrostatic equilibrium is attained. It is these converged hy- 
drostatic models that form the basis of the parametrization 
described below, even though the hydrodynamical calcula- 
tion evolves away from this initial static structure towards a 
steady fiow solution. We however check retrospectively that 
our parametrization also applies to the evolved structure 
(see Section 5). 

Under the assumption that the thermal and hydrody- 
namical timescales are decouplecjf] a coupling of the MO- 
CASSIN code with a hydrodynamical code is trivial in princi- 
ple. It simply consists of feeding the density fields produced 
by the hydrodynamical code into MOCASSIN to calculate gas 
temperatures that are then in turn fed back to the hydrody- 
namical code. While conceptually easy, this type of coupling 
is computationally prohibitive. 

In order to render the calculation tractable, we exploit 

The thermal timescale in the flow is typically less than a year 
and thus more than an order of magnitude less than the hydro- 
dynamic time scale. 
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the fact that the wind will be optically thin to high energy 
photons. In the density and temperature regimes expected 
in a photoevaporative wind, radiative heating and ioniza- 
tion compete with collision-initiated cooling and recombina- 
tion and hence the ionization equilibrium stro ngly depends 
on th e ionization parameter ^ = Lx /nr'^ (see iTarter et alJ 
1 19691 ). The ionization parameter can be conceptually con- 
sidered as the ionizing energy available to a gas particle at 
a distance r from the source and will determine the ioniza- 
tion state of the gas and hence th e gas temperature at that 
point (e.g. llgea fc GlassgoldlllOg^ ). We show in Figure [2]the 
mean relation between ionization parameter and tempera- 
tu re in the X-ray heated region of the hydrostatic models 
of lErcolano et ahl (|2009lj ) (upper panel), where we define 
the X-ray heated region as comprising gas for which the 
absorbing column to the central star is less than 10^^ cm~^ 
(roughly the penetration depth of 1 keV photons). The scat- 
ter of the MOCASSIN grid cells about the relation shown is 
small (~ 0.3 dex) thus showing that the ionization param- 
eter is indeed a good determinant of the gas temperature. 
This confirms the central importance of X-ray heating in 
determining the thermodynamics of the heated atmosphere 
(note that a tight relationship between ionization parame- 
ter and temperature is not expected in general in the case of 
heating dominated by EUV radiation). In Section 4 we will 
consider the case of photoevaporation from discs containing 
inner holes for which we need to improve our parametriza- 
tion in the lower temperature regime. The lower panel of 
Figure [2] depicts the result of a similar calculation of a disc 
which is truncated at an inner radius of 16 AU; we have 
confirmed that this provides also a good fit to models with 
different inner hole sizes and hence use this relation in all 
the hydrodynamic models in Section 4. Note the similarity 
of the relations for the identical ionization parameter range. 
This implies that it is the direct field that is dominating in 
both cases. This should be contrasted with the pure EUV 
case, where the diffuse field dominates for the optically thick 
primordial disc and the direct field dominates for discs with 
inner holes. 

Since these relations show that the gas temperature 
varies monotonically with ionization parameter in the Xray 
heated regime, we can thus use the ionization parameter 
(which depends only on the incident flux and local num- 
ber density) in order to prescribe the gas temperature for 
regions lying at columns of < 10^^ cm"'^ from the central 
source. 

Fo r columns larger than lO^'^cm"^ lErcolano et al.l 
(|2009bl 'l (ECD09) found that the gas and dust were ther- 
mally coupled. Hence we use the dust temperatures of 
iD'Alessio et all Hool) to fix the gas temperatures for 
columns greater than lO^'^cm"^. In order to check that our 
temperature parametrization is appropriate for the final con- 
verged flow solution, we feed back this density field into MO- 
CASSIN and compare the temperatures thus obtained; these 
results are presented in ^S] 

2.2 Radiation-Hydrodynamic Methods 

In order to accurately determine the mass-loss rates and the 
kinematic and morphological structure of the evaporative 
flow, we need to compute the temperature distribution in 
the flow and feed this into a hydrodynamic code. We have 
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Figure 2. Plot of Temper ature as a function of io nization param- 
eter for FS0H2Lxl model llErcolano et al.ll2009bl) with numerical 
parametrization shown in red. The black asterisks show the mean 
temperature in the ionization parameter bins obtained from MO- 
CASSIN, with a scatter about the mean of ~ 0.3 dex. 

modified the publicly ava il able zeus-2d hydr odynamic code 
(|Stone fc Norman|[l992al lbl: IStone et al.lll992 l to include the 
heating from X-ray and EUV radiation, by writing a new 
module that calculates the gas temperature using the meth- 
ods discussed in i|2.1l We then use this temperature to up- 
date the internal energy density of gas (u) in the wind via 
the ideal equation of state: 

u = ?^pr (1) 

where ks is the Boltzmann constant, is the hydrogen 
mass and 7 is the ratio of the specific heats. The mean molec- 
ular weight is fixed in the flow, to a value of /i = 1.37125. 
The internal energy is only required to update the pressure 
needed for evaluation of the momentum equation. The pres- 
sure is thus determined again from the ideal gas law, via: 

p = (7 - \)u (2) 

and hence the procedure is independent of choice of 7. This 
thermal update is done during both the source step and 
transport step in order to minimise the time to a converged 
solution. 
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We assume azimuthal and mid-plane symmetry and 
since we expect a photoevaporative flow to be approxi- 
mately spherical at large radius, we use a spherical grid. We 
neglect magnetic fields and self- gravity. The rotation op- 
tion in ZEUS-2D is enabled, which introduces the necessary 
pseudo-forces arising from rotating frames. The inner and 
outer radial boundaries are set to outflow, while the angu- 
lar boundaries are given the appropriate symmetry bound- 
ary conditions. We choose to use the standard van-Leer 
(second-order) interpolation scheme and the von-Neumann 
& Richtmyer parametrization for the artificial viscosity with 
gvisc = 2.0. 



3 PRIMORDIAL DISC 

In order to determine the flow structure for the disc for the 
majority of its standard viscously accreting phase we run a 
single hydrodynamic model. We expect the structure of the 
flow to be set by the temperature of the X-ray heated re- 
gion and, to a lesser extent, by the temperature of the disc 
atmosphere which is controlled by optical heating of dust 
in the disc's surface layers. The structure of the underlying 
disc is insensitive to this provided that it is very optically 
thick, both in the X-ray and optical regions. This is amply 
fulfilled by primordial (untruncated) discs over many orders 
of magnitude in accretion rate, thus justifying the fact that 
we use a single hydrodynamic model to evaluate the wind 
mass loss profile for the entire period of the disc's evolution 
prior to gap opening (see Sections 4 and 6). A similar ar- 
gument has previously been applied also to the EUV case, 
where the mass loss profile is unaffected by the total col- 
umn in the disc, provid ed it is optically thick in the EUV 
l|Hollenbach et al.lll994D . 
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Figure 3. The density structure of a primordial disc; with the 
initial disc shown in (a) and the 'steady-state' flow solution in 
(b). The red line indicates the r=2/3 surface from the central 
star in both plots. In (b); streamlines are plotted at 5% intervals 
of the integrated mass loss rate, with the sonic surface plotted 
in magenta. Note the smoothness of the sonic surface indicates 
that we are extremely close to a true 'steady-state'. The snapshot 
shown is at a simulation time of 10 orbital periods at lOOAU. 



3.1 Numerical Method 

The grid is chosen to be in the range 6 — [0,7r/2] and 
r = [0.33, lOO]At0 with 100 grid cells in the angular di- 
rection and 250 in the radial. Grid cells are logarithmically 
spaced in the radial direction, so that we have adequate res- 
olution at small radius: in particular we ensure that we can 
resolve the onset of any fiow and can resolve the scale height 
of the X-ray 'dark' disc at all radii. As initial conditions, the 
vertically hy drostatic density and temperature structures 
generated bv lErcola no ct al. (20093) were interpolated onto 
the ZEUS-2D grids and assigned Keplerian rotation (suitably 
adjusted for the small radial pressure gradient). The disc 
structure was then allowed to evolve hydrodynamically using 
the modified code, with the gas temperatures being updated 
as described in i|2.1l at every timestep, until a 'steady-stat^fl' 
disc/wind system was obtained. It should be noted that the 
density structure of the X-ray 'dark' regions is also allowed 
to evolve hydrodynamically, i.e. it is not reset to the original 
density structure at every time-step (compared with previ- 
ous EUV hydrodynamic simulations e.g. I Alexander! l2008bl : 



^ Note throughout this paper we use (r, 9, <j>) and {R, ip, z) to 
distinguish between spherical and cylindrical polar co-ordinates 
respectfully. 

* Where the term 'steady-state' strictly means quasi-steady-state 
or more formally M^, X t -C M^isc ■ 



iFont et al]l2004l which - since they did not model the EUV 
'dark' disc - reset the base density of the flow to its original 
value at each time-step). We emphasise that the mass loss 
over the duration of the simulation is much less than the 
initial disc mass. 

The evaporative flow was found to have converged after 
~4-6 orbital timescales at the outer boundary (convergence 
was also checked by determining the Jacobi potential was 
conserved along streamlines: See fjS]). The flow structure af- 
ter 10 orbital timescales at the outer boundary is used in 
our analysis. 

3.2 Results 

In Figure |3] we show the density structure of the wind and 
underlying primordial disc with streamlines overlain. The 
Figure shows that the flow originates near the point at which 
the X-ray flux can penetrate, i.e. a column of lO^'^cm"'^. 
This also indicates that the flow launching surface is signifl- 
catly higher than the disc's photosphere as seen by the star 



(red solid line i n 3(b) where we have used the opacities of 
iD'Alessio et al.ll2001l to determine this surface). Also, com- 
paring the intial disc structure (Figure |3(a)[ | and the final 
flow solution (Figure 3(b) I it is clear that the optically thick 
disc region is very similar and is thus observationally con- 
sistent with the SED of classical T-tauri stars by construc- 
tion, with the surface of optical depth to stellar photons of 
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Figure 4. Cumulative mass loss rate as a function of cylindrical 
radius, with the numerical fit used in the viscous evolution shown 
as the solid line. Note the smooth profile with only one 'bump', 
indicates the fiow is extremely close to a true 'steady-state'. 



T = 2/3 lying at ~ 0.1-R. We stress that although the X-ray 
heated wind launches from well above the disc mid-plane, 
both the wind and the upper regions of the underlying disc 
are optically thin to stellar photons and thus the fraction 
of young stars that are optically obscured b y their discs re- 
mains at the ~ 10% level implied by the lo'Alessio et al.l 
lj200ll ) models. 

Figure |3] also shows that the mass loss is dominated 
from the region 5-20 AU, although there is still significant 
mass-loss at larger radii. We note that the streamlines are 
close to radial at large radius similar t o those obtained 
from pure EUV models l|Font et al.l |2004| : [Alexander et al.l 
l2006al : lAlexandedbOOSbl ). The observation that the flow is 
approximately radial at large radius allows us to be confi- 
dent that the outer boundary has not affected our results, 
sinc e radial supersonic outflo w is treated exactly in ZEUS- 
2d IStone fc NormanI Il992al ). The total mass loss rate is 
found by integrating the mass flux across a spherical sur- 
face situated at a radius of 85% of the outer grid boundary 
(sufficiently far from the outer radial boundary to avoid any 
effects from spurious reflections due to non-radial outflow); 
and X-ray 'dark' gas is explicitly excluded from our measure- 
ment. We then determine the cylindrical radius from within 
which the integrated mass loss rate occurs by following the 
streamlines down to the base of the flow. Figure!?] shows the 
cumulative mass-loss rate as a function of cylindrical radius, 
and shows that the mass-loss rate is dominated by the region 
between 5 and 20 AU, but significant mass loss occurs also 
at larger radii with 50% of the mass loss occuring between 
18-70AU. This should be compared to the pure EUV case 
where only ^10% of the total mass-loss occurs outside 18AU 
l|Font et al.ll200j ). The X-rays cause a much broader wind 
profile by being able to heat gas at larger radii than pure 
EUV irradiation due to their greater penetrative power. The 
total integrated mass-loss rate over the entire disc was found 
to be 1.4xlO"*M0yr"^ 



Table 1. Table listing the input parameters for simulations of 
inner hole discs, where i?i„ is the inner edge of the X-ray 'dark' 
disc 



Model 


-Rin (AU) 


Grid range (AU) 


A 


8.3 


[1.33,90] 


B 


9.7 


[1.67,90] 


C 


11.7 


[3.35,100] 


D 


14.2 


[3.35,100] 


E 


17.7 


[3.35,100] 


F 


21.1 


[3.35,100] 


G 


30.5 


[3.35,100] 


H 


64.5 


[10,450] 


I 


86.9 


[15,500] 



4 INNER HOLE FLOWS 

As will be shown in i)6.2l photoevaporation combined with 
viscous evolution results in a gap opening in the primordial 
disc; in contrast to the case of planet clearing, where streams 
cross the gap at certain azimuths, in the present (axisym- 
metric) case the inner disc is starved of re-supply once the 
gap is opened. The inner disc then drains onto the star leav- 
ing the outer disc open to direct irradiation by the central 
star. The fiow expected from such a system is of course dif- 
ferent from a primordial disc's photoevaporative fiow, and 
in this section we set out the methods used to calculated 
mass-loss rates from discs with inner holes. 

4.1 Numerical Method 

Given we expect the mass-loss rates to vary with inner hole 
size as seen in the EUV case l| Alexander et al]|2006al ). we 
have run several models with inner holes at different radii, 
out to a radius beyond which disc holes become hard to 
detect observationally (i.e. at around ~ 100 AU). We used 
a regularly space spherical grid in order to maximise reso- 
lution while minimising computational effort. The grid was 
constructed with 400 cells in the radial direction and 200 
cells in the angular direction, again sufficient to resolve the 
onset of any fiow and the scale height of the 'dark' disc at all 
radii. All other zeus-2d parameters were kept the same as 
above. Table [T] summarises the input parameters used. As 
initial condition, the hydrodynamic primordial disc model 
was cut at the inner hole radius. Then the density profile was 
smoothed over several pressure scale lengths at inner edge, 
the we adjusted the angular velocity to match the new im- 
posed radial pressure gradient. The model was then allowed 
to evolve hydrodynamically until it attained a steady fiow 
solution. 

4.2 Analytic Model 

Before we perform detailed hydrodynamic simulations to 
determine the fiow structure of transition discs it is use- 
ful to develop a theoretical framework in which to analyze 
the simulations. In order to estimate the mass-loss r a te scal - 
ing we adopt a similar method to [Alexander et al.l (|2006ah 
and assume the mass loss rate per unit area is oc pvi{R) 
where vi{R) is the local launch ve locity of the ioni z ed gas . 
However unlike the calculations of [Alexander et all (|2006ah 
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Figure 5. The 'steady-state' density structure of a transition 
disc and its photoevaporative wind with an inner hole of 14.2 
AU (model D), the streamlines are plotted at 5% intervals of the 
integrated mass loss rate, with the sonic surface plotted in red. 
Note the sonic surface being close to smooth indicates that we 
are extremly close to a true 'steady-state'. The snapshot shown 
is at a simulation time of 10 orbital periods at lOOAU. 



where the gas was assumed to be isothermal, we have no 
such restriction in our calculations. We expect streamlines 
that are launched from the disc's inner rim (which will turn 
out to dominate the mass loss) to start on a nearly radial in- 
ward trajectory (see Figure lXSjl . Inevitably such a trajectory 
will result in a centrifugally driven wind (whose force term 
scales as l/R"^) resulting in a launch velocity that scales as 
the escape velocity Ve{R). Furthermore, if we assume that 
the effective emitting area of the wind scales as the square 
of the inner disc radius (Rin) then we can write the wind 
mass-loss rate as scaling with pve{R)Ri„ oc ppi^l^ . We will 
use the code to determine the number density at the base 
of the flow at the disc's inner rim (ni„) and check out our 
estimate that in this case 



M 3/2 



(3) 



The actual scaling of M with however needs to 
be determined numerically since the ionization parameter 
at the flow base (and hence ni„) depends on the tempera- 
ture at the inner rim and this needs to be calculated using 
the radiation-hydrodynamic algorithm for each hole size. We 
parametrize the dependence of mass loss rate on i?i„ in Sec- 
tion g^] 



4.3 Results 

All inner hole flows share the same general morphology: the 
mass-loss is dominated by regions close to the inner edge, 
but significant mass-loss still occurs at larger radii. Figure 
O shows the flow structure of model D with an inner hole at 
14.2AU. We extract mass-loss rates in an identical fashion 
to the primordial disc case and note again that the radial 
structure of the streamlines at large radii allows us to be 
confident that the outer boundary has not affected our mass- 
loss rates. 
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Figure 6. Plot showing analytic model {Mvjiud I i^ia ^f/^) 
agrees well with the simulations. Error bars shown indicate er- 
rors in riin obtained by averaging over cells at the inner edge. 
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Figure 7. Plot showing M at different inner hole sizes; The blue 
point represent mass-loss rates from models A-I, while the black 
line shows the numerical fit used in the viscous evolution. 



Figure [6] shows the ratio M/ni„ reproduces the pre- 
dicted R^^^ scaling discussed at the end of >14.2I i.e. it is 
consistent with our assumption that the effective launch area 
scales as _R?„ and that the launch velocity scales as the es- 
cape velocity. 

In Figure [7] we present the mass-loss rates as a func- 
tion of inner hole radius complete with the numerical fit 
used in the viscous evolution. This shows that for _Ri„ out to 
~ 30 — 40AU the mass- loss rates fall approximately as R~,^'^ 
while at larger i?i„ the mass-loss rates begin to increase. We 
can understand this behavior by noting that equation ^ 
and the definition of the ionization parameter implies that 
Af oc so that the scaling of photoevaporation rate 

with radius depends on the radial variation of f and hence on 
the form of the temperature versus ionization parameter re- 



lation. From Figure 2(b) we see that this contains a point of 
inflection at log^ — —5.2, logT = 3.2 which is encountered 
when the Ri„ is around 40 AU. At smaller radii and hence 
higher temperatures the variation in ionization parameter is 
small for a given temperature variation and hence the mod- 
est decline in the flow temperature with inner hole radius is 
achieved at nearly constant ^. Beyond Rin ~ 40 AU, a mod- 
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est decline in temperature corresponds to a larger decline in 
^ and hence the fall off in M becomes less steep and even 
starts to increase somewhat with increasing 7?i„. It should 
however be stressed that the radial variation of M with _Ri„ 
is very mild (see Figure (O). 



5 TESTS 

In order to be confident that the solutions obtained are ac- 
curate and are unaffected by the assumptions discussed in 
i)2.2l we have conducted various tests which we discuss be- 
low: 

(i) A numerical convergence test using double the resolu- 
tion confirms the flow is converged to ~5%. 

(ii) In order to be confident the parametrization of X-ray 
heating is correctly reproducing the temperature structure 
in the fiow we have re-calculated the gas temperatures of 
the steady-state solution using mocassin. Figure |S] shows 
a comparison between the temperatures determined by the 
ZEUS-2D algorithm and those obtained using MOCASSIN. This 
shows that we have errors <15% for 50% of all grid cells 
and errors <40% for 95% of the grid. In Figure |9] we show 
the comparison between the gas temperatures determined 
using the zeus-2d algorithm and those obtained using MO- 
CASSIN at various cuts in cylindrical radii. It shows that 
we have good agreement throughout the grid especially in 
the launch region of the flow (region of large temperature 
gradient) and in the regions of the flow which are subsonic 
(red regions) where the thermal pressure gradient can be a 
first order effect. We note that the disagreement between 
the temperatures in the flow/disc transition region is an ar- 
tifact of MOCASSIN which imposes a hard transition between 
the molecular and atomic zones (see discussion in ECD09), 
hence the errors quoted can be considered an overestimate. 

(iii) It has been assumed that the gas and dust are ther- 
modynamically coupled in regions that are optically thick 
to X-ray irradiation. In order to test this assumption we 
have performed a run where we allow the gas in the X-ray 
dark regions to evolve adiabatically. We then compare the 
timescale for the gas an d dust to return to thermal equilib- 
rium through collisions ijWhitworth fc Clarkd IiQQtI ) to the 
local dynamical timescale and we flnd that the dynamical 
timescale is ~ 100 the thermal timescale at the minimum. 
Thus we can conclude that the gas and dust will be ther- 
modynamically coupled when not penetrated by the X-ray 
flux. 

(iv) In the following sections, we assume that the mass- 
loss profile obtained for the primordial disc is valid for the 
majority of the discs lifetime (i.e. the viscously accreting 
phase) . During this period the disc structure will change as 
material drains onto the central star. In order to check our 
mass-loss rate is is valid for a varying (optically thick) disc 
structure, we compare it to the mass-loss rates obtained for 
the adiabatic runs mention above, which produce markedly 
different (adiabatic) disc structures and find good agreement 
between the models. 

(v) In our hydrodynamic modeling we have used the dust 
temperatures (and h ence gas temperature o f the X-ray 'dark' 
region) obtained bv iD'Alessio et al] (|200ll ). Since the dust 
temperatures will be linked to the disc structure and vice- 
versa, and since the disc structure has evolved from the orig- 



inal hydrostatic equilibrium structure, is im portant to assess 
how a ccurately the dust temperatures of the lD'Alessio etahl 
(l200ll ) model can describe the dust temperatures in the hy- 
drodynamic solution. We've calculated the dust tempera- 
tures of both the original model and the flow solution using 
MOCASSIN and find we have agreement to within the Monte- 
Carlo errors. 

(vi) In order to determine whether our models are close 
to a true steady-state we compute the Jacobj^ potential 
along the streamlines, noting that this will be constant in a 
steady flow. For each streamline we determine the enthalpy 
term by paramatrizing the run of pressure and density as 
a baratropic relation (i.e. P{p)) along the streamline, then 
integrating along the streamline. However we note that, the 
flow is not strictly baratropic or irrotational and thus don't 
expect the Jacobi potential to be constant for all stream- 
lines. Figure [To] shows that the Jacobi constant is conserved 
along the streamlines. We flnd that the Jacobi constant is 
conserved to within 3-4% indicating we are very close to a 
true steady-state flow. 

(vii) ECD09 models do not include molecular cooling in 
the thermal balance. While (as shown in Figure 3(b) I the 
flow is launched in a region higher than the optically thick 
disc, the lack of molecular cooling could in principle be a 
potential source of error on the gas temperature in the disc 
interiors, and as such the caus e of a hydrostatically 'p uffed 
up' disc. Although as shown by Glass gold et al.l (|2004l . their 
Figure 4) molecular cooling only comes in at large column 
densities and, even there, it only has a very limited role 
when compared to dust-gas collisions, it is useful to test 
whether a cooler (flatter) disc interior would result in signif- 
icant changes in the photoevaporative mass loss rate of the 
system (due to neglecting molecular cooling, or some other, 
as yet unidentifled cooling mechanism). We can do this by 
considering the mass-loss rates obtained from the adiabatic 
models described above (point iii & iv). The use of an adia- 
batic equation of state results in a cooler disc structure and 
hence a lower launching height. While the exact quantitative 
result is sensitive to disc heating, the overall mass-loss rates 
in all models were found to be ~ 10"* Mq yr"^ This effect 
is small since the change in gravitational energy between 
the midplane and the launching height is of order 10% of 
the total gravitational energy gained by the wind, as shown 
in Figure [14] 



6 VISCOUS EVOLUTION 

The evolution and dispersal of an irradiated protoplanetary 
disc is controlled by the simultaneous actions of the photoe- 
vaporative wind and viscous torques in the disc. A numerical 
simulation of disc dispersal must therefore simultaneously 
account for photoevaporation and viscous evolution. 



6.1 Method 

The diffusion equation for the evolution of the sur- 
face density is (|Lvnden-BeU fc Pringljll974l : jPringld [l98ll : 



^ The Jacobi potential is the analogue to the Bernoulli potential 
in a rotating frame. 
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Figure 8. Histogram of cells with temperature errors less than a 
given value. Black solid line refers to errors in cells in entire grid, 
red dashed line to errors from cells that are locally sub-sonic (in 
the rotating frame). 



Figure 10. Plot of the Jacobi potential along various streamlines 
in the range 2-50AU, plotted against the normalised streamline 
length. The variation in the mean value along each streamline is 
~ 3% 
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Figure 9. Radial cuts showing temperature comparison for 
R=2,10,20AU. Solid line indicates temperatures calculated by the 
ZEUS-2D algorithm, where the red solid line indicates where those 
cells are locally subsonic (in the rotating frame) and in the pho- 
toevaporative flow. Blue stars represent the temperatures deter- 
mined by MOCASSIN. 



[Alexander et al ] |2006bl ): 



as 
dt 



RdR 



(4) 



Given that t h e in itial disc structure taken from 
iD'Alessio et al.l (Il998l ) has a surface density that scales 
approximately as E oc R~^ , an initial surface density 
structure of the form 



E(i?,0) 



Md(0) 
2-kRiR 



exp{—R/Ri 



(5) 



is used, where Md{0) is the initial disc mass and Ri is 
the initial scale size of the disc. This initial structure 
is a similarity solution of t he viscous diffusion equation 
(jLvnden-BeU fc Pringljll974l ) when the viscosity scales lin- 
early with radius and can be parametrized as 



and the mass-loss in the wind can be neglected. The vis- 
cous scaling constant vi can be determined from (see 
'Hart mann et all Il998l : IClarke et al] l200ll : [Alexander et all 
2006b|)l 



M{t = 0) 



3MdiO)ui 
2Rl 



(7) 



V = uiR/Ri 



(6) 



The integration of equation Q is performed using a first- 
order explicit scheme, the grid consists of 1000 points equi- 
spaced in R'^ and spans the rang e [0.0025, 2500] AU with 
zero-torque boundary conditions (|Pringle et al.lll986l ). The 
radius beyond which of the disc mass is initially lo- 
cated (> Ri) is taken to be 10 AU and the initial accretion 
rate is 5 x 10~^Mq yr"'^. While we adopt a viscosity scaling 
oc i? in order to be consistent with the R~^ surfa c e den sity 
profile of our i nitial disc structure, IClarke et al.l ()200ll ) & 
Alexander et al . (2006b) considered the effect of using a dif- 
ferent viscosity parametrization and concluded that, while 
the absolute timings would change, the overall qualitative 
nature of a photoevaporative 'switch' was unchanged. 

The relative timescales between the three stages (vis- 
cously evolving, gap opening and disc clearing) are indepen- 
dent of the initial disc mass (for a fixed initial accretion rate) 
and can thus be scaled accordingly for any chosen value. In 
order to set a reasonable absolute value for the disc lifetime, 
we have chosen an initial disc mass of O.2M0, which (given 
the initial accretion rate) corresponds to an initial viscous 
timescale at the disc's outer edge of 2 x lO^yr; this in turn 
corresponds, for our chosen case of Ri = 10 A.U. to a viscous 
alpha parameter of ~ 10~^ at Ri. 



6.1.1 Surface Mass-loss rates 

The functional forms of Eu, used are those obtained for the 
primordial disc models and inner hole models. Since we ex- 
pect the inner hole to drain on a short timescale, as a first 
approximation we either use the primordial disc wind pro- 
file, or the inner hole mass-loss profile and do not attempt 
to model the transition in any detail. The transition oc- 
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Figure 11. Snapshots of the surface density plotted at 
0%, 20%, 40%, 60%, 70%, 75%, 77%, 79%... 99% of the disc's 
total lifetime (see footnote 4), The surface density is in dimen- 
sionless form and in units of A'/^(0)/27riJ^ (see equation [St. 



curs once the radial column density to the inner edge is 
< 10^^ cm-^ 

6.2 Results 

Figure [TT] shows the evolution of the surface density of the 
disc. At early times the evolution of the disc is similar to a 
standard viscous accretion disc with a power-law decline in 
surface density and accretion rate. Once the accretion rate 
is comparable to the mass-loss rate due to the wind (which 
occurs at ~ 60% of the disc lifetim^fl), the photoevaporative 
wind begins to affect the surface density profile in the range 
5-30 AU; a gap then opens at lAU after 78% of the disc life- 
time. Once the gap opens, the inner disc quickly drains on a 
timescale of order the viscous timescale at 1 AU. During this 
draining time the inner edge of the outer disc is continually 
photoevaporated outward due to the significant mass-loss at 
larger radius. 

Once the inner disc has drained, the outer disc is ex- 
posed to the direct photoevaporative wind, and the inner 
edge of the disc has now moved to 9AU. The total mass 
remaining in the disc after the inner disc has drained onto 
the star was 9% of the intial disc mass and the outer disc 
is then dispersed entirely in the remaining 15% of the disc's 
lifetime. We find, in total 84% of the intial disc mass is ac- 
creted onto the central star, with the remaining being lost 
in the wind. 

In Figure [12] we show the accretion rate onto the star 
throughout the disc's lifetime and show that the disc evolves 
much faster compared to a standard viscously accreting disc. 
We also see that, even prior to gap opening, the disc spends 

^ At this point we define the disc lifetime as being the total life- 
time when photoevaporative clearing is taken into account; i.e. 
the time to clear out to 87AU (beyond which M is not accurately 
known); Due to the power law nature of viscous draining, it is 
hard to define a 'lifetime' of a control disc (with no photoevap- 
oration). However, to be indicative: the model disc that we have 
run has a total lifetime when photoevaporation is included of 
« 1.3 Myr; if such a disc were simply viscously evolving without 
photoevaporation it would become optically thin in the region of 
the disc that dominates the K band emission after 0(100) Myr. 
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Figure 12. Plot of accretion rate onto the star against time. 
The black solid line shows the model for the X-EUV wind, while 
the blue dashed line sh ows the pure EUV canonical model (see 
[Alexander et al.ll20063 . Figure 1). 

a significant portion of its lifetime (~ 15%) with the accre- 
tion rate onto the star being less than the photoevaporative 
wind rate: the accretion rate onto the star is lxlO~^M0 
yr"'^ at the point that the gap opens, which is about an 
order of magnitude less than the wind mass loss rate. This 
effect (whi ch has been terme d 'photoevaporation starved ac- 
cretion ' bv lDrake et"al]|2009t ) is a consequence of the radially 
extended mass loss in X-ray heated discs so that, prior to 
gap opening, there is a significant region of the disc where 
the accretion rate decreases towards small radii. This behav- 
ior was not seen to be a prenounced effect in previous pure 
EUV models since, in the EUV case, the mass loss is more 
strongly localised and the gap opens in a disc that is hardly 
depleted by the wind at large radii. 



7 COMPARISON WITH PREVIOUS WORK 

Several other authors have considered photoe vaporative flow 
driven by either, EUV or X-ray radiation and lErcolano et aU 
(2009b) have considered the combined effect of X-EUV ra- 
diation. The predicted mass loss rates due to pure EUV 
photoevaporati on have been acc urately determined as 1.1 x 
lO"^°M0yr"^ (|Font et al.ll2004l) for a O.7M0 star and an 
ionizing flux of 10*^ s"^ Also [Alexander et all (|2006a') found 
that when the inner hole first forms, the mass loss rate rises 
by about a factor ten and then increases mildly as the hole 
grows according to B^i'^ ■ 

As discussed in Section 1, previous estimates of X- 
ray photoevaporation rates to be found in the literature 
are qualitatively different from the values presented in the 
present w ork since they are not b ased on hydrodynamic cal- 
culations. [Alexander et al.l l|2004) estimated mass-loss rates 
per unit area at ~ 20 AU which were similar to those 
obtained in this work. E CD09 argued for similar surface 
mass-loss rates as those of [Alexander et al] (|2004[ ). but es- 
timated a much larger launching regions (i.e. over a re- 
gion 5-50 AU). They therefore concluded that X-ray pho- 
toevaporation would dominate over EUV photoevaporation 
at all but the lowest of X-r ay luminosities. The estimates of 
[Gorti fc HollenbachI l[2009l ) suggested much lower mass loss 
rates than the EUV case, implying that X-rays would only 
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have an indirect effect on disc dispersal. ECD09 include a 
detailed discussion of the differences in the thermochemical 
calculations of Gorti and HoUenbach 2009, which may have 
led to their result, and concluded that the rather harder X- 
ray spectrum employed by Gorti and HoUenbach may have 
been an important contributory factor. 

The mass loss rate of 1.4x lO~*M0yr~^ that we obtain 
from our radiation hydrodynamical calculations can be com- 
pared to the estimate of ECD09 for the same input param- 
eters (model FS0H2Lxl of ECD09) of 6.7x10"'' Mgyr'^ 
This difference arises from the fact that ECD09 assumed 
that the gas escaped at the sound speed at the point at 
which the internal energy of the gas was equivalent to the 
escape potential. In reality the flow is launched sub-sonically 
from deeper in the disc where the density is higher. It is nev- 
ertheless reassuring that the two values are of a similar order 
of magnitude. We also flnd that the mass-loss occurs over a 
large radial range of 1-70AU, in contrast to the EUV case 
where it is steeply peaked about one characteristic radius 
(for a O.7M0 star this would be 6.2AU). Given the mass 
loss rates obtained and the monotonic temperature : ioniza- 
tion parameter relations found (see Figure [S]), it is clear in 
this case that X-ray photoevaporation will be a dominant 
effect in disc evolution. 

Coupling the mass-loss profiles to secular viscous evo- 
lutions models has only been previously attempted for 
the pure EUV case by [Alexander et al.l l|20 06b) where 
the mass-loss profile had been accurately determined. 
In order to compare ou r mode l with what the EUV 
model of [Alexander et al.l (l2006bh would predict with the 
same input parameters , we take the canonical model (see 
[Alexander et al.l l2006b'. Figure 1) and multiply both the ini- 
tial disc mass and the initial disc viscous timescale by a fac- 
tor 4, thus ensuring that the rescaled model has the same 
initial accretion rate, disc mass and viscous timescale as the 
present model. □ In the (rescaled) EUV case, it was found 
that the disc dispersed entirely after ~24Myr with the tran- 
sition stage occurring on a short ~ lO^yr timescale (approx- 
imately 3% of the total disc lifetime). Our results are con- 
sistent with this 'switch' nature of disc dispersal (albeit in 
somewhat softened form: see Figure [12]). In the present (X- 
EUV) models we find that a gap in the disc opens at ~1AU 
after 78% of the discs lifetime (i.e. nominally at an age of 
~lMyr). Once the inner edge of the outer disc is exposed to 
direct radiation we find the outer disc is quickly removed on 
a timescale ~ lO^yr with a mass loss rate scaling as -R"'' '^. 

While X-ray driven photoevaporation produces the 
same qualitative results when it comes to disc dispersal, 
there are several distinct differences. Prior to gap opening 



' Ignorance of the magnitude of the viscosity in discs means that 
the absolute values for disc lifetimes that we quote here in the 
EUV and X-ray heated case should not be used as grounds for 
favoring one model or another: we have scaled the model used 
here so that its total lifetime is an observationally reasonable 
~lMyr and the EUV lifetime is then rather long, but we could 
have just as easily scaled down the initial mass and initial viscous 
timescale in such a way that the EUV lifetime was observationally 
reasonable and the X-ray lifetime was unreasonably short. The 
purpose of our present comparison is therefore to look at the ratio 
of timescales in X-ray and EUV models with the same underlying 
viscous parameters. 
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Figure 13. In this Figure we plot the surface mass-loss rate for 
the X-EUV wind (solid black line) obtained in this work and th e 
EUV wind (dashed blue line) llFont et al.l2004l : lAlexanderll2008bl) . 
We have scaled the peaks of the profile to be the same value in 
order to compare the relative broadness of the different profiles. 
Clearly this shows the X-EUV profile is much broader than the 
EUV profile. 



the surface density begins to fiatten, indicating that the ac- 
cretion rate is dropping below that of a steady disc (i.e. that 
the disc is in a state of 'photoevaporation starved accretion' 
as discussed above). This effect is much less obvious in the 
EUV case owing to the more restricted radial extent of the 
EUV wind. Once the gap opens at around a few AU (simi- 
larly in both cases), the inner edge of the outer disc evolves 
rapidly to ~ lOAU in the X-ray case due to the fact that this 
region was already significantly depleted, due to the broad 
mass-loss profile. In the EUV model, the inner edge remains 
at lAU during disc draining, due to the steep mass-loss pro- 
file and then quickly moves outwards to ~ 20AU as the pho- 
toevaporation rate rises by an order of magnitude when the 
inner edge of the disc is directly irradiated. The difference 
in behavior of the edge of outer disc during inner disc drain- 
ing, can be readily understood when the surface mass-loss 
profiles are compared for the EUV case and X-EUV case. 
Figure [12] clearly shows the X-EUV profile is much broader 
than the EUV profile; the EUV wind thus has very little ef- 
fect on the disc beyond Vg whereas the X-EUV wind subjects 
the disc to a period of depleted accretion inwards of 25AU 
prior to gap opening. During the clearing of the outer disc, 
the mass-loss rate falls as i?"" "^ (out to 30-40AU) in the X- 

1/2 

ray case, compared to increasing as i?i„' in the EUV case. 
This difference is due to the fact that in the EUV case, where 
the ionized region is approximately isothermal, the fiow is 
launched at roughly constant (~ sonic) velocity and a sim- 
ple Stroemgren sphere argument implies that the density in 
the launch zone scales as _R~^'^; in the X-ray case, where a 
range of temperatures are present, the effective launch speed 
from the exposed inner rim of the disc is the local escape 

— 1/2 

speed which scales as R^^ . The density in the launch zone 
(which corresponds to a nearly constant value of the ioniza- 
tion parameter) falls off slightly less steeply than Both 
these effects mean that EUV photoevaporation is relatively 
favored at large radius. We estimate that the radius at which 
the direct EUV photoevaporation rate becomes comparable 
to the X-ray photoevaporation rate is > 200AU and hence 
the EUV will remain unimportant in the case of discs around 
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solar type stars with X-ray luminosities '-^lO^^erg and 
greater. It may however be relevant for systems with lower 
X-ray luminosities. Given the large scatter of X-ray lumi- 
nosities for stars of similar masses, a full discussion of the 
parameter space is beyond the scope of this paper, but will 
be considered in future work. 



8 MODEL LIMITATIONS 

While this model is an important step forward in un- 
derstanding disc dispersal, it is important to discuss the 
limitations of our models. We have neglected the ef- 
fects of non-ionizing FUV radiation in our calculations. 
iGorti fc HoUenbachI (|2009l ) have shown that FUV radiation 
is capable of heating gas to above the escape temperature at 
large radius (100-200AU), and predict mass loss rates of a 
similar magnitude to those obtained from our X-EUV mod- 
els. Given that these mass-loss rates are far more uncertain 
than those obtained from a self-consistent hydrodynamic 
model we cannot make any definite predictions of the ef- 
fect FUV inclusion would have on our model. Qualitatively, 
however, the existence of a further FUV wind at large ra- 
dius would starve the inner disc of material and would thus 
accelerate the final inside-out clearing described here. 

Throughout our modeling we assume that the gas is al- 
ways in thermal equilibrium with the radiation based on the 
comparisons between the thermal and flow timescales dis- 
cussed in ^2.11 This assumption can be tested by comparing 
the gained energy flux of the wind to the ionizing luminos- 
ity. Figure [14] shows a plot of specific energy for the primor- 
dial disc. Considering the streamlines around the median 
mass-loss value (launching from 15-20AU) the total change 
in specific energy along these streamlines from launch (with 
a negative specific energy of ~ 10^^ erg g"'^) to exit from 
the grid (with a positive specific energy of ~ 5 x lO^'^ erg 
g~^) is ~ 6 X 10^^ erg g^^. Given the grid boundary is far 



from the sonic surface (see Figure 3(b) I this can be consid 



ered close to the maximum gai ned energy flux of the wind 
(c.f. Parker wind: |Parkeilll958l '). In order to compare with 
the input energy rate (Ltot), we can compute the mechanical 
luminosity of the wind using L,nech ~ x Ae (where e is 
the speciflc energy) which evaluates to ~ 10^^ erg s~^. Thus 
it is clear that the energy gained in the wind is much less 
than the input energy from the X-EUV spectrum (4 x 10^° 
erg s~^). For all flows (both the primordial and inner-hole 
flows) we explicitly calculate the gained energy flux in the 
wind and flnd that its <^ Lx (~<8%). Hence the assump- 
tion of thermal equilibrium is valid in this case. Hence the 
terms in the energy equation from advection and pdV work 
can be neglected. 

During the transition between a primordial disc and a 
transition disc, we assume the mass-loss proflle for the for- 
mer, until the column to the inner edge of the outer disc is 
less than lO^^cm"'^. In reality we would expect the transi- 
tion to be a smooth process. Given the difference in mass-loss 
rates between the primordial and transition disc models is 
small, and the transition time is rapid, such an approxima- 
tion will have little effect on the overall evolution. However 
a more detailed approach will be required when it comes to 
predicting the statistical distri bution of inner hole size s as a 
function of stellar parameters (jAlexander et al.ll2006bh . 




Radius (AU) 



Figure 14. Figure showing the specific energy of gas particles 
throughout the grid. Streamlines are plotted at 5% of the to- 
tal mass-loss rate. The dark blue lane indicates the transition 
between a overall negative specific energy (bound gas) and an 
overall positive specific energy (unbound gas). This clearly shows 
that the total energy gained by gas launched into the wind is 
less than Lx with a total efficiency ~ 8% compared to the input 
X-ray energy. 



9 OBSERVATIONAL CONSEQUENCES 

We have shown that the X-EUV photoevaporation model de- 
scribed here can reproduce the 'two-timescale' evolution ob- 
served in young solar type stars, in agreement with current 
observational constraints on relative timescales. In this sense 
it is qualitatively similar to the EUV switch model as set out 
in [Alexander et al] ('2006b'): in particular the second (rapid 
clearing) timescale is similar and so is the radius (~ 1 AU) 
at which the gap first opens. The temperature distribution 
in the flow is however rather different in the two cases (the 
X-EUV wind is somewhat cooler: compare Figure]^ with the 
EUV case where the wind is close to 10*K throughout). It 
is thus an immediate task to compare the predictions of the 
new model with observational diagnostics that have been 
studied in the EUV case (Ercolano et al in prep.). In partic- 
ular, the strength and l i ne profile of the [Ne II] emission at 
12.8;U i n (|Pascuccill2007l: iPascucci et al.1 '2009': 'H erczeg et al.1 



I2OO7I : iNaiita et al.l l2009l : iFlaccomio et al.. ,200gf ) has been 
sug gested as an observational test for X-ray irradiated discs 



e.g.l'GlassEold et al. 1120071 : iMeiierliik" et al ]|2008l : lGudeill2008l : 

[Schisano ct al. 200y); our calculations open up the possi- 
bility o f detailed modeling, analogous to that of I Alexander] 
()2008bl ') for the EUV only case. Likewise, the X-EUV model 
needs to be assessed with regard to the production of for- 
bidden line emission (including 0[I]) as in the study of 
Font et al. (2004) for the EUV case. 

The important difference in terms of disc secular evo- 
lution between the X-EUV and the EUV only model is the 
magnitude of the photoevaporative mass loss rate which is 
larger in the X-EUV case by 2 orders of magnitude. The 
wind mass loss rate is now comparable with t ypical accre- 
tion rates measured in Classical T Tauri stars l]Natta et al.] 
]2006f ). This means that whereas the EUV switch provided a 
mechanism for the final shutdown of a protostellar accretion 
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disc when its accretion rate had fallen to very low (nearly un- 
detectable) values, the X-EUV wind (at least for the rather 
high X-ray luminosity employed here) should cut off a disc 
in its prime. This raises a number of issues: 



9.1 Accretion rates in T Tauri stars and their 
relationship with X-ray luminosity 

If all discs were subject to the vigorous X-EUV winds found 
in this model, then young stars should rarely exhibit accre- 
tion rates less than a few times 10~^Mq yr~^, since sub- 
sequent draining is rapid (see Figure I12p . This is clearly 
untrue observationally, as this accretion rate is close to the 
median for Classical T Tauri stars. However it needs to be 
borne in mind that T Tauri stars exhibit a ~ 2 orders 
of m agnitude scatter in in X-ray luminosities at a given 
mass jPreibisch et al.ll2005h and that the value that we have 
adopted here is modestly above the median value. We have 
not yet undertaken full hydrodynamical calculations with a 
range of X-ray luminosities but the hydrostatic estimates of 
ECD09 tentatively suggest that the photoevaporation rate 
may scale roughly linearly with X-ray luminosity (cf the 
EUV case which, being recombination limited, predicts a 
mass loss rate that scales as the square ro ot of the EUV 
output of the star; iHollenbach et all ll994l V Therefore we 
would expect that in stars with a lower X-ray luminosity 
the plot of accretion rate versus time (equivalent to Figure 
I12p would steepen at a lower accretion rate. 

What is the expected sign of the dependence of accre- 
tion rate upon X-ray luminosity within the subset of stars 
that have not yet cleared their discs (i.e. in Classical T Tauri 
stars)? Some care is needed in phrasing this question. If we 
have a coeval cohort of stars, then evidently the stars that 
are stronger X-ray emitters would be expected to exhibit 
lower accretion rates, as a consequence of 'photoevapora- 
tion starved accretion'. On the other hand, if we have a 
mixture of (disc bearing) stars of all ages then in a statisti- 
cal sense the strong X-ray emitters will be younger (because 
the predicted disc lifetime is less) and they will have higher 
accretion rates on average (because the low luminosity X-ray 
sources live for a longer time and spend most of this time 
with accretion rates modestly above the (low) wind rate). 

A clear observational consequence of the lower disc 
lifetime in luminous X-ray sources is that such sources 
spend a higher fraction of the pre-main sequence state 
in a disc-less (Weak Line T Tauri) state: consequently 
one expects that Weak Line T Tauri stars should be 
more luminous X-ray sources on average than their disc 
bearing (Cla ssical T Tauri) counterparts. This is indeed 
the case (e.g. Neuhauser et al.|[l995l : iFlaccomio et al.ll2003l : 
iPreibisch et al.l 20051 ) though this has also been variously 
explained in ter ms either of X -ray absorption in accretion 
columns (Gregorv et al.l |2007| ) or confine ment of the X- 
ray p roducing corona in accreting systems l|Preibisch et al.l 
|200a) (in other words, it is conventionally assumed that ac- 
cretion suppresses X-ray emission rather than the effect pre- 
dict ed here whereby X- ray emission suppresses accretion: see 
also iDrake et al.l |2009| for a discussion of this possibility) . 
Naturally, this effect cannot be quantified until we have ver- 
ified through hydrodynamic simulations how the wind mass 
loss rate scales with X-ray luminosity. 



9.2 Disc lifetimes 

The fact that the X-EUV photoevaporation rate is much 
higher than that for EUV only has obvious implications for 
disc lifetimes. For example, we find that for the same initial 
parameters, the disc lifetime is around 1.5 orders of magni- 
tude lower in the present case than when compare d with the 
EUV only calculations of lAlexander et al] l|2006bl ). However 
it should be borne in mind that in both cases the absolute 
lifetime depends on the the time for viscous evolution down 
to a critical accretion rate and this therefore depends both 
on the initial disc mass and the viscosity in the disc. The 
magnitude of the latter can in any case be adjusted so as 
to reproduce any reasonable disc lifetime. What will be in- 
teresting, however, is to examine the distribution of relative 
disc lifetimes that results for a plausible dis tribution of disc 
initia l properties and X-ray luminosities (c.f. lArmitage et al.l 
l2003h . 



9.3 Transition discs (inner hole sources) 

We have shown that X-ray luminous T Tauri stars should 
evolve through the following sequence of states: a) optically 
thick (untruncated) stage dominated by viscous evolution b) 
inner gap phase, with a low density inner disc and accretion 
on to the central star that is less than the wind mass loss rate 
( < 1O~*M0 yr~^) and c) inner hole phase with an evacu- 
ated inner disc and essentially no accretion onto the star and 
an outer disc mass of < O.OIM©. This is qualitatively sim i- 
lar to the EUV case described in [Alexander et al] (|2006bl ): 
important quantitative differences are the (~ 100 x) higher 
upper limit on the accretion rate in stage b) in the X-EUV 
case and also the somewhat higher upper limit (roughly an 
order of magnitude) on the outer disc mass during stage c). 
It has often been argued that the relatively high accretion 
rates and outer disc masses in observed inner hole sources is 
a counter-a rgument against the relevance of (EUV) photoe- 
vapo ration (jNaiita et al.ll2007l : rAlexandeij|2008al : lKim et all 
I2OO9I ). We thus see that X-EUV photoevaporation alleviates 
both of these problems. 

It would however be unwarranted to necessarily con- 
clude that all (or even the majority) of observed T Tauri 
stars with spectral evidence for cleared inner holes are cre- 
ated in this way. Indeed it should be stressed that even 
if we judge that the bulk of such sources are created by 
other mechanisms, this does not diminish the important 
role of photoevaporative clearing in the majority of discs. 
It needs to be recalled that the evolutionary sequence de- 
scribed above implies - in the absence of other interven- 
ing processes (such as the formation of planets or grain 
growth) - that all disc systems should evolve through a 
short-lived inner hole phase. The observed frequency of in- 
ner hole sources (of order 10% of all disc-bearing sources) 
is broadly compatible either with a scenario in which all 
observed inner hole sources are a short-lived evolutionary 
phase undergone by all stars or one in which the majority 
of observed inner hole sources are instead a long-lived evo- 
lutionary phase undergone by a subset of all T Tauri stars. 
There is no shortage of other mechanisms proposed for in- 
ner hole creation, inclu ding for example clearing by planets 
or binary companion s (iQuillen et al.l |2004| : IVarniere et al.l 
I2OO6I : iRice et al.|[20o3 . |2006| ). grain~growth in the inner disc 
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Dullemond fc DominikI I2OO5I : ICieslal |2007| ). photophoresis 
Krauss et al.l 2007h and inside out evacuation by the mag- 



netoroational i nstability ([Chiang fc Murrav-Clavll2007l') fsee 
e.g. reviews bv lAlexandejboOSal and discussion in lKim et al.l 
I2OO9I ): what is notable is that these alternative mechanisms 
are not linked to a mechanism for clearing the residual outer 
disc and are thus not associated with a short-lived phase 
prior to complete disc dispersal. Therefore although they 
may apply to some of the observed inner hole sources they 
cannot represent an evolutionary norm for young stars with- 
out over-predicting the incidence of inner hole sources. 

Given that we cannot distinguish on statistical grounds 
between a short-lived transistion phase experienced by all 
stars or a long-lived phase undergone by a minority, we 
need instead to look at the properties of individual sources, 
bearing in mind that the properties of these sources will 
partly depend on the observational criteria used for object 
selection. Therefore, for examp l e, the sample of inner hole 
sources compiled by iKim et all (|2009l ) on the basis of a dip 
in the Spitzer IRS spectrum (in the range 5 — 40/im) are 
largely a group of objects with high accretion rates on to the 
central star 1Q~^Mq yr~^) and high outer disc masses 
O.OIM0). Whereas the X-EUV models are close to being 
able to reproduce these conditions in a short period briefly 
after gap opening, they should not be representative of pho- 
toevaporative holes in general. In particular, X-EUV models 
are not compatible with large (> 20 A.U.) holes that are 
vigoro usly accreting (such as SZ Cha, GM Aurig a or UX 
TauA ; iNaiita et a"Lll2007l. ICalvet et alll2005. . Espaillat et all 
20071). On the othe r hand, the inner hole sources compiled by 



Cieza et al.l (|2008l ) are qualitatively different, being associ- 



ated with lower disc masses and low or undetectable levels of 
accretion: these may indeed include sources with photoevap- 
orating holes in the evoluti onary stage c) above. (Like wise 
the mm survey of Andrews I Andrews fc Williams! [ioOSi con- 
tains objects that may be in this category, although none 
of these have the high ~ O.OIM© disc masses that are pre- 
dicted by X-EUV models at the onset of outer disc clearing) . 
It should be stressed that EUV photoevaporation models 
are mildly biased towards larger hole sizes and that this 
effect is somewhat stronger in X-EUV models, due to the 
larger outer disc masses and the gentle decline in wind rate 
at large radii (Figure [13]). We would therefore expect that 
photoevaporating holes would be more strongly represented 
amongst sources with large inner holes as one would expect 
to uncover with Herschel or ALMA. 

A final point about the formation mechanism for ob- 
served inner hole sources is that it is a robust prediction 
of photoevaporation models that there should be no flow 
across the gap once it is opened. Thus following gap open- 
ing, the draining inner disc is not re-supplied and the ratio 
of its mass to the accretion rate onto the star should be 
a good measure of its lifetime in the draining phase. Ob- 
served inner hole sources are generally characterised by the 
dust content of their inner discs and in most cases the gas 
mass is unknown. However, me asurements of CO emission 
in TW Hydra and GM Auriga feal vk et aLll2007h indicate 
very low inner disc gas masses, which, when combined with 
the measureed accretion rates in these objects would imply 
improbably short lifetimes in the absence of re-supply. This 
argument therefore is a strong indication that the inner disc 



is re-supplied in these objects, as is a natural expectation in 
planet clearing models. 



10 CONCLUSIONS 

In this work we have presented the first self-consistent hy- 
drodynamic model of an X-EUV driven photoevaporation 
wind. The mass-loss rate of 1.4x lO~*M0yr~^ occurs over 
a large radial extent of 1-70AU. This should be contrasted 
with the much lower ~ IQ^^'^ Mq yr^^ values obtained from 
EUV photoevaporation alone and with the fact that in the 
EUV case the mass loss rate is more concentrated at a scale 
of a few AU. 

When combined with viscous evolution, X-EUV photo- 
evaporation produces a qualitatively similar pattern of inner 
hole creation followed by inside out clearing as is found in 
EUV models. The initial size of the inner hole (about an 
AU) and the rapid clearing following gap opening is simi- 
lar in X-EUV and EUV models. The important difference is 
that this occurs earlier, and at much higher accretion rate, 
in the X-EUV case. This implies that inner holes can be 
associated with more massive outer discs and with higher 
accretion rates on to the star than was possible in the EUV 
case. Our results thus weaken the objections made previ- 
ously to EUV photoevaporation with regard to inner hole 
creation. For example, the X-EUV model can open a gap 
even in the case of massive outer discs such as that in TW 
Hydra. 

Nevertheless there are generic predictions of photoevap- 
oration models (whether EUV or X-EUV) that appear not 
to be compatible with features of some of the best stud- 
ied inner hole sources. Namely: i) photoevaporation models 
predict that gas and dust in the inner hole, as well as di- 
agnostics of accretion onto the star, should be observable 
only during inner disc clearing and that during this phase 
the ratio of inner disc gas mass to accretion rate should be 
of order the lifetime of this phase. In TW Hydra and GM 
Auriga, however, the observed ratio is very low and would 
imply an unacceptably short lifetime of this phase. This low 
value probably implies that, in contrast to what is predicted 
in photoevaporation models, the gas and dust in the inner 
hole is continuously replenished from the outer disc, ii) Pho- 
toevaporation models also predict that the phase of inner 
hole draining is associated with small holes (i.e. of < 10 AU 
scale). There is no mechanism to produce large holes that 
contain gas and undergo vigorous accretion onto the star (as 
in GM Auriga). 

On the other hand, both EUV and X-EUV models can 
explain the detection of outer discs in around 10% of stars 
which show no evidence for inner disc emission/accretion 
diagnostics (| Andrews fc WiUiamj [ioosh . Future measure- 
ments with Herschel or ALMA are likely to uncover a class 
of objects with large inner holes and we would expect this 
population to be dominated by photoevaporative holes. 

Finally, we stress that the X-EUV case investigated here 
corresponds to the case of a rather large X-ray luminosity 
(2 X 10'"' erg s~^) and further work is required in order to ex- 
plore the effect on disc evolution of the wide range of X-ray 
luminosities measured in T Tauri stars. A clear qualitative 
prediction is that disc lifetimes are expected to be shorter 
in X-ray luminous sources and that these sources should be 
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be over-represented among disc-less (Weak Line) T Tauri 
stars. This prediction is in line with observations: though 
the larger X-ray luminosity in Weak Line T Tauri stars is 
conventionally ascribed to supp ression of X-r ays in accret- 
ing systems, we here follow Dra ke et alj (|2009l ) in suggesting 
that another factor is that X-rays may instead suppress ac- 
cretion. 
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